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Abstract 

This work is devoted to the study of integration with respect to 
binomial measures. We develop interpolatory quadrature rules and 
study their properties. Local error estimates for these rules are derived 
in a general framework. 
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1 Introduction 

In this paper we develop quadrature formulae for the numerical integration 
with respect to (w.r.t.) binomial measures. The binomial measure fi a , where 

< a < 1 is a parameter, is a probability measure on an interval of the real 
line, say [a, b], that is characterized by the following (self-similar) property: 
let / a dyadic subinterval of [a, b] and bisect I in the left and right parts 

1 = I L U I R ; then 

Moi(/r) = Oi/jL a (I) . (1) 

When a = 1/2 we trivially obtain the probability measure proportional 
to the Lebesgue measure on [a, b]. Without loss of generality we will fix 
[a, b] = [0,1]. The family of measures {n a }a has important features that 
makes interesting the study of their properties. Firs of all, it is a fam- 
ily of pairwise mutually singular (see definition 6.7 [23J) Borel measures: 
fjL ai -L fJ-a 2 if a i a 2 ; thus in particular each \x a is singular w.r.t. the 
Lebeasgue measure Va ^ 0.5. Moreover, \x a is a continuous measure, i.e. 

MM) = o Vxg [o,i]. 

Binomial measures naturally appears in problems related to to the proba- 
bility theory of sequences of independent trials, and in particular to the so 



L i±l 

2 fe' 2 k 



C [0, 1[, where j,k G N. 



called Bernoulli proces^]. Let Ik ■ 
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email {calabro corbo} @unicas.it 

Each independent trial has two possible outcomes -success or fail- with fixed proba- 
bilities, respectively a and (1 — a). 
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Then fj, a (Ik) is exactly the probability that in the first k trials of the 
Bernoulli process we have a number of successes equal to the number of 
Is in the binary expansion^! of j. 

Also, following the notations in [15], the measure \i a can be defined as the 
unique measure that satisfied the following balancing equation: 

J fdfi a = (1-a) J f(x/2) d^ a + aj f(x/2 + 1/2) d/z« (2) 

V/ 6 C 1 . This relation is usefull when multifractal properties are studied, 
see [TTJ [22] for an introduction and [3[ [H [191 HO] for possible developments. 



The construction of quadrature rules for integration w.r.t. binomial mea- 
sures is considered for two reasons. 

The first is, obviously, for the calculation of the integrals because for these, 
if a 7^ 1/2, we cannot describe the solutions using an analogous of the funda- 
mental theorem of the integral calculus. Notice, moreover, that the relation 

gives also that the functional L[f] = J f dfi a is a "refinable linear func- 
tional" , as recently defined in [17] , with the Bernoulli shifts as stretch-shift 
operators and mask [2(1 — a), 2a]. Thus the calculation of these integrals 
seams to be an active research problem. 

The other reason is that when we exhibit a quadrature rule we are giving 
in an implicit manner a way, via combination of Dirac delta, for a decom- 
position of the measure. Ad example, if the moments are preservecQ, the 
inverse problem is known as moment problem and has been considered also 
for balanced measures, see [H EJ Q2]. Notice also that the use, for a gen- 
eral measure, of a decomposition that involves the calculation of quadrature 
rules based on balanced measures has been explored recently, see [18] . 

In this paper we analyze how polynomials can be integrated on dyadic 
intervals and introduce interpolation based integration rules. 

The work is organized as follows. In the second section we report some 
analytical background and write some useful technical identities. In the third 
section we introduce quadrature with respect to the measures \i a . Finally, 
we list some remarks and possible future work. 



2 Preliminary Results 

In this section we present a self-contained introduction of the binomial mea- 
sures. For the sake of completeness we also list some results that we will 

2 A detailed study of Hausdorff dimension of sets related to the averages of binary digits 
is performed in [3 OH]. 

^Notice that the functions (x/2) and (x/2 + 1/2) are usually called Bernoulli shifts. 

4 This feature is, as we will see, easily related to the degree of exactness of the quadrature 
rule. 
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Figure 1: The normalization of the density function 2 k f ajk with a = 0.3 
and k = 6 



often apply in the sequel. First of all we recall the mean value theorem. 

Proposition 2.1 (Mean value theorem). Let /j, be a positive finite measure 
defined in [c,d], and consider (f> G C°([c, d]). Then there 3£ G [c,d]: 

1 f d 

MM]) Jc 

For sake of clearness, and without loosing of generality, we will restrict 
ourself to the case of the interval of integration to be / = [0, 1]. An important 

role will be played by dyadic intervals X k = 

j < 2 k . We will call k the order of the interval. 

Let Ha t k be the probability measure with constant density on dyadic intervals 
of order k given by: 

f a , k (x) = 2 k a n ^(l-a) k ~ n ^, A< x< ^±i ) j = 0...2 fc -l 

where: 

n(j) = #{ls of the binary expansion of j} . 

Two useful relations involving this measure are summarized in the fol- 
lowing lemma. 

Lemma 2.2. Given h, k £ N, the following hold true: 

fl a , k (Xj) = Va,h( X j) k>h 

^, k+h {X k +\.) = ^ k {X k )^ h {X*) Vi = . . . 2 h - l,j = . . . 2 k - 1 . 



L i±l 

2 fc ' 2 k 



kN, < 
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Proof: The first relation is immediate, while the second is easily proved 
noticing that, due to the fact that since i runs from . . . 2 h — 1 : 

n(j2 h + i) = n(j2 h ) + n(i) = n(j) + n(i) => 

a n{j2h+i) {l - a )k+h-n(j2 h +i) = [ a n(i)(! _ a )*-»(j)] [ a «W (1 _ a ) h ~ n ^} . 

□ 

We can now prove the convergence for fixed a of the sequence of measures 

{^a,k}k- 

Proposition 2.3. The sequence of measures {fJ> a ,k}k converge in the weak- 
star sense (see definition 1.58 in JBj). 

Proof: Since the set of probability measures is weak-star compact, we 
only need to prove existence of the limit of the sequence: 

1 i>(x)dfi aik ) V0ec°([o,i]) . 
J k 

We will prove that this is a Cauchy sequence. Given e > 0, since (j) is 
uniformly continuous, there exists hy G N such that: 

sup \4>{x) - 4>{y)\ < e . 

\x-y\<2- k 

Now, take k,2,ki > ko and consider: 



x) dfj, ajkl - / 4>{x) dfi a ,k 2 



< 



2 fc 0-l 



3=0 



2 R 0-1 



X k 



x) dfj, aM - / 4>{x) d/i a! k 2 



x k 



3=0 

2 fe 0-l 2 fc 0-l 

= E w^) - ^,k (4°) < E (4°) = e . 

j=0 i=o 

where we have applied the mean value theorem (£ij,£2,j are points in X* ) 
and the first statement in lemma 

The sequence, thus, converges pointwise, as requested. □ 

Definition 2.4 (Binomial measures fj, a ). Fixed a G (0, 1), we will call bi- 
nomial measure /i a the weak-star limit measure of the sequence {n a k}k- 

We will denote, as usual, with L v ^ a the space of the p— integrable func- 
tions with respect to [i a . 

We state now a "change of variable" type result for measures fi a . 
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Lemma 2.5. Let f be in L^. Then, for each dyadic interval Xf we have 



that: 



f 

JO 



f(x) dfi c 



^--/ f(2 k x-j)dfi c 
X 1 ") Jx* 



Proof: Let us first consider the case where the integrand / is the char- 
acteristic function \E of a measurable set E C [0,1]. By lemma [2^21 the 
formula is true if E is a dyadic interval, and by summation for a finite union 
of such sets. Taking the supremum of such kind of functions we can obtain 
the formula for characteristic function of open sets; moreover taking once 
again the infimum we can get the formula for any xe and by linear combi- 
nation for any simple function. 

Let us now take / 6 L , / > (in the general case we can write / = f + — f~ 
where / + and /~ are respectively the positive and negative parts). Since 
/ is the pointwise limit of a monotonic sequence of simple functions (see 
theorem 1.17 on 1 231) we obtain the result. □ 



Note that lemma 12.51 cannot be extended to any affine change of variables, 
i.e. in general ft f{x) dfi a + Ma( ^ a+Aj) f* +X /(^) dfi a . 

Now let us see how to calculate the moments. Take s G N, we apply lemma 
I to write: 



x s d\i c 



1 



o 



1 



a 



1/2 



(2x) s dn c 



a 



1/2 



x s dfx a 



(3) 



Now, we can notice that: 
-i 



x s d\x a 



o 



1/2 



x s d\x a + I x s dfi 

1/2 



applying another time Lemma 

»l/2 



f 1 / 2 a f 1 ' 2 
/ x s d[i a + / (x + l/2) s dfi. 

Jo 1 - « Jo 



1/2 



x s dfi a + 



a 



1 



a 



1/2 



E 

q=0 



r s-q_ 



a 



1/2 



X S dfi a + 



a 



1 



a 



1/2 



E 

9=1 



2" 



dfi c 



2<-i 



a 



1/2 



X S dfi a + 



n 



1 



E 

9=1 



1 

21 



1/2 



— / X 



s-q 



dfl a 



dfi a 



5 



applying the relation ([3j) in the parenthesis: 



- I s rl 

" 1 - a J X 1 ~ a q Ti L 2" 2^ J 

Summarizing, we have that: 



s\ 1 1 — a 



x s q d/j, c 



g=l 



(4) 



Now, substituting in the second term of (J1J) the relation seen in ([3]) we have 
that: 



x s d\x c 



1 I- a 



1 — a 



and therefore we get the following: 

Proposition 2.6. Moments of the measures fi a are connected by the fol 
lowing recursive relation: 



x s dfi a 



9/ Jn 



□ 



We will use piecewise s-polynomial interpolation in the next chapters. 
For this reason we explicitly calculate the integrals of polynomials in the 
dyadic intervals Xj applying lemma [2~51 



3 Quadrature rules 

In this section we want to introduce the numerical integration rules. We will 
call integration rule a choice of p + 1 distinct points ( q S [0, 1] (called nodes) 
and of values (5 q (called weights). Let f(x) G L^, we will denote b)H: 

q=0 

5 Note that, beside the integrability condition, we will always apply quadrature rules 
to functions with a finite number of discontinuities and everywhere defined. 
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Table 1: Integration Rules 



rule nodes ( g weights (3 q 



Wi 



a 



degree of 
exactness 



a 1 1 

8a + 3 V-264a 2 + 264a + 81 \ 18" ~ 9 

14 42 2 2V-264a 2 + 264a + 81 3 

8a + 3 V-264a 2 + 264a + 81 1 , 18a - 9 

14 + 42 2 2V-264a 2 + 264a + 81 

/ a(l — a) 



1/2 2 if a ^1/2 

/a(l-a) 1/2 3 if a = 1/2 

a + - : 



Notice that, as pointed out in the introduction, this is equivalent to consider 
as an approximation of the measure ji a the following combination of Dirac 
delta: \i a « Ei?=oA? ' S C 9 - 

We will call degree of exactness of the formula I p with respect to \x a the 
greatest positive integer r such that the considered decomposition maintains 
the same moments up to order r: 

[ x q dfx a - I p (x q ) = \/q < r , q G iV • (5) 
Jo 

This relation in the case of r = gives a re- normalization on the weights: 

<7=0 

If we consider to replace the function with a polynomial, we can calculate 
the integral in an exact manner by means of the seen formulae for monomi- 
als, proposition 12,61 If the polynomial is chosen as the one interpolating the 
function in the nodes, this leads to the so called interpolation-based inte- 
gration rules. Fixed the nodes, the weights of the formula can be calculated 
integrating, as in the case of the Lebesgue measure, the so called Lagrange 
fundamental polynomials (see equation (9.2) in [21]); in our case this inte- 
gration is to be made with respect to the measure \x a - 
A well known theorem, valid for general positive measures (see [H]), states 
that this rules can give degree of exactness up to r = 2p + 1. It is also well 
known that, fixed the nodes £ 9 , q = . . .p there exists an unique choice of 
weights that leads to a formula of degree r > p, and this rule is necessarily 
interpolation-based. 

Optimal rules (from the degree of exactness point of view) are usually 
referred to as Gauss-Christoffel, see [13] for a complete review in the general 
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case of a positive measure. We have calculated the first two of them with 
respect to the measure \x a and these are are summarized in table [TJ denoted 
by G p . It can be seen how already G% is quite laborious to describe [f| for 
these measures. 

If we consider p+1 equi-spaced points the corresponding interpolation-based 
quadrature rules that we obtain are called Newton-Cotes formulae. The first 
of these rules are summarized as NC P in table [21 We can notice that the 
weights are negative for some choices of a at the end-points. Ad example, in 
the case of NC2 the weights are positive in the case of a £ (1/4,3/4), while 
if a < 1/4 the last weight is negative and in the case of a > 3/4 the first 
is negative. Moreover, these formulae are of degree of exactness^] p, that is, 
obviously, the lowest for interpolation-based rules. 

We can moreover fix some parameters and use nonlinear relations ([5]) to con- 
struct the highest degree interpolatory formula satisfying these constraints. 
As example, we can construct a formula with two equal weights and of degree 
of exactness 2, see rule Wi in tabled) This formula, although, is meaningful 
only if a £ [1/4,3/4] (otherwise nodes are outside of the integration interval) 
and coincides with Gi in the case a = 1/2 (for shortness in the sequel we 
will refer to the case a = 1/2 as the Lebesgue case). 

Another useful example of three point interpolation-based rule is in table [5J 
denoted by GL 2 . This formula is obtained fixing the two extreme nodes and 
considering the only rule of degree of exactness 3 with three nodes. This 
procedure leads to the so called Gauss Lobatto formulae. Notice that for 
this formula the weights are always positive and that in the Lebesgue case 
GL 2 = NC 2 . 

3.1 A-priori error estimates 

In this section we develop a-priori error estimates for interpolation-based 
rules. The first that we will see relays on the corresponding formulae for 
interpolation errors. Recall that if Il^f (a?) is the (unique) p-polynomial 
interpolating function f(x) at the nodes Q q , q = 0, . . . ,p and / is sufficiently 
regular, we have that (see [2T] equation 8.7): 

f(x) - n /jfg (x) = — - * u p (x) for some £ x £ [0, 1] , 

6 For the construction of these rules, the classic procedure is to calculate the nodes as 
zeros of a polynomial of degree p+1 determinate by means of an orthogonal procedure 
(weights can be calculated consequently). For a survey about this or other procedures, 
for example involving the so called Jacobi matrix, we refer to [131 114j . 

7 We have p+1 only in the case a = 1/2 for even p. 
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Table 2: Integration Rules 2 



rule 

NC 
NCi 

NC 2 
NC, 



NC 4 



GL 2 



nodes ( q 

1 

2 



1 


1/2 

1 


1/3 
2/3 

1 


1/4 
1/2 
3/4 

1 



(3a + 2)/7 
1 



weights /3„ 



1 — a 

a 

4a 2 - 7a + 3 
-8a 2 + 8a 
4a 2 - a 



/ 



/ 



1 

315 



9a 3 + 24a 2 - 22a + 7 \ 
27a 3 - 51a 2 + 24a 
-27a 3 + 30a 2 - 3a 
\ 9a 3 - 3a 2 + a 
( 256a 4 - 992a 3 + 1572a 2 - 1151a + 315 \ 
-32a(32a 3 - 94a 2 + 99a - 37) 
"■ 128a 2 + 73a -9) 



a: 



\ 



24a (64 

-32a(32a 3 - 34a 2 + 9a - 7) 
a(256a 3 - 32a 2 + 132a - 41) 
I (a-l)(5a-6) ^ 

(3a + 2) 
98a(a - 1) 



/ 



(3a + 2)(3a-5) 
a(5a + 1) 

\ (5 - 3a) J 



degree of 
exactness 

if a + 1/2 

1 if a = 1/2 

1 

2 if a ^ 1/2 

3 if a = 1/2 



4 if 1/2 

5 if a = 1/2 



where u> p is the so called nodal polynomial: 



g=0 



Applying this estimate we obtain: 

f 1 f(x) dfi a - i p (f) = f 1 [f(x) - n /)C , (x)] dfi a 

Jo Jo 

= ^h)\Sj {p+1) ^ x)Up[x)dlla 



(6) 



In the particular case that / is a polynomial of degree p+k with 1 < k < p+1 
then / — Rf£ q can be factorized as u p (x)q(x) where q(x) is a polynomial of 
degree k—1. Therefore if u p is L^ a orthogonal to all polynomials of degree 
up to k — 1 we obtain that the rule has degree of exactness p + k, that is a 
well known result for quadrature with respect to a general positive measure, 
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see theorem 2.1 in [14 
Let us write a simple formula for the error estimate involving both the 
derivatives of order p+1 and p + 2. Call k>p~(x) = max{0, u p (x)}, uj~(x) = 
max{0, -ujp(x)}, K+ = J^uj^(x)d^ a and K~ = f u>~ (x) d\i a then by © 
we get: 



f(x) dfl a - Ip(f) = 


' f f {p+1) (^H( x ) dua - C f^ +1 \^H(x) dn a 
Jo Jo ' . 

/^ +1 )(£i) f u+(x)d^ a -f(P +1 \^) f\-(x)dfi a 



1 

1 

(P + 1)! 
1 

1 

(P + 1)! 



^C[/ (p+1) (£i) - / (p+1) (6)] - K - k+] / (p+1) (e 2 



^/ (p+2) (C 3 )(Ci-6) + / (p+1) (6) / w P (x)dM« 



4H) 



In ([BED we have applied proposition 12,11 simply noticing that /^ p+1 H^x) = 

n p (a:)] can |_ )e regaled as a continuous function on [0, 1]. 
Note that in the same way we could put in evidence K~ instead of equa- 
tion (|6E1) obtaining an analogous estimate. This estimate is quite difficult 
to use, because the K a can be explicitly computed only in few cases (one of 
these cases can be seen in the remark [3JJD . 



Another useful estimate relays on the Taylor expansion of /. Consider I 
be of degree of exactness r and consider / to be r + 1 times derivable. Take 
the Taylor expansion of the function up to the power r in the point x = 0: 



fix) = P5(x) + R r +1 (x) where: 



r+l. 



n=0 



(r + l)! 



.r+l 



I 



Take now the error of the formula I: 

J f(x)dfi a -I(f) 

P^(x) dfl a - I(P r ) + J R r Q +1 (x) dfia ~ HR ' 
f (r+l) {0 



r+l\ 



(r + l)! 



„r+l 



dpL c 



N 
q=0 



r+l 

g 



(7) 



Notice that we have applied a discrete version of the mean value theorem as 
in (|6b| . In this estimate an important role plays the error committed when 
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calculating the first moment where the formula is not exact. For this reason, 
in tables [HUMS] are reported errors made by some of the constructed formulae 
when calculating the moments of the measure. In particular, in each table 
is fixed a. Errors are quite significant, as it can be seen comparing with the 
exact values reported in the first line of the tables. 

We test some of our formulae also in the case of a polynomial of forth grade 
and on the 20-th moment: 



5x 4 + 6x 3 



h 



..20 



(8) 



In figure [5] we plot the exact errors made from the four formulae NC2, NC3, 
GL2 and Gi when a varies in (0,0.5]. 

Remark 3.1 (Formulae of degree p+1). Notice that we have constructed, in 
the the cases of odd nodes, two formulae, namely Go and GL2 that generalize 
the corresponding equispaced rules and maintain the degree of exactness p+1 
that is achieved in the Lebesgue case. In the five point case it seems to us 
that there is no straight generalization of NC4 having degree of exactness 5 
for all a. In particular, if we fix the choice of the external points to be as in 
the equispaced case and move the midpoint in function of a, the rule exists 
only for some choices of a (in an interval of 1/2 approximately of amplitude 
0.3185846387J, and the new point is: 

320a 3 + 512a 2 - 237a - 158 
31(64a 2 - 64a + 9) ' 

It is also possible to construct a formulae on 5 points of degree of exactness 
5 moving the second and forth point in a rigid linear manner with respect 
to a, but we obtain a rule that does not coincide with NC4 in the Lebesgue 
case. The nodes are: 




7^/31 11 . 
I I 1/2 

27900 465 30' 1 ' 



1939 7\/31 11 
27900 + 465 + 30' 



We will see in RemarkW. 
5 that coincides with NC4 
cases, six points. 



how to construct a formula of degree of exactness 
in the Lebesgue case but that uses, in the other 

□ 



We conclude this section with two examples on how to apply these esti- 
mates (|6|)- (l6Hp . the three point Newton-Cotes rule NC2 (where the estimate 
is quite pessimistic) and the formula GL2 (where we obtain a good estimate). 
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Figure 2: Plot of the true errors (Jg 1 / dfj, a — over the a. Test function 
are in equation (JSj) . 

Remark 3.2 (Error Estimate 1: NC2). In the case of the rule with 3 equi- 
spaced points we can calculate the constants in (IMl) simply noticing that: 



k: 



k: 



uj 2 (x)>0 if i£ [0,1/2] 

u 2 (x) = -u 2 (x - 1/2) ifxe [1/2,1] 

1/2 , , , . , . . . , a(l-a) 2 (4-a) 

a; 2 (x) d> Q = (1 - a) / w 2 (x/2) d^Q, = — 
Jo 



28 



1/2 



w 2 (a;) d^Q, 



a 



x / 2 , , , a 2 (l-a)(4-a) 
w 2 (x) d\i a = — 



1 - a jo 
a(l - a) (4 - a)(l - 2a) 



28 



28 



Taking the absolute value we obtain: 



f(x)dfi a -nc 2 f(x) < 

|a(l-a)(4-a)(l-2a)| 

168 u 1 

ail — a) (4 — a) min(a, 1 — a) . , 

H max r 

168 ^ 
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Table 3: Calculated moments, case a = 0.05 





s= 1 


2 


3 


4 


5 


X s du . 5 


0.0500 


0.0183 


0.00860714285714 


0.00448142857143 


0.00248679339478 


NC (x s ) 


0.5000 


0.2500 


0.12500000000000 


0.06250000000000 


0.03125000000000 


Nd(x s ) 


0.0500 


0.0500 


0.05000000000000 


0.05000000000000 


0.05000000000000 


NC 2 (x s ) 


0.0500 


0.0183 


0.00250000000000 


-0.00541666666667 


-0.00937500000000 


EC 3 (x s ) 


0.0500 


0.0183 


0.00860714285714 


0.00591798941799 


0.00539021164021 


<S a (x s ) 


0.0500 


0.0025 


0.00012500000000 


0.00000625000000 


0.00000031250000 


Gi(x s ) 


0.0500 


0.0183 


0.00860714285714 


0.00407148526077 


0.00192634556203 



Table 4: Calculated moments, case a = 0.3 





s=l 


2 


3 


4 


5 


X s du . 3 


0.3000 


0.1600 


0.10200000000000 


0.07136000000000 


0.05300129032258 


NC (x s ) 


0.5000 


0.2500 


0.12500000000000 


0.06250000000000 


0.03125000000000 


NCi(x s ) 


0.3000 


0.3000 


0.30000000000000 


0.30000000000000 


0.30000000000000 


NC 2 (x s ) 


0.3000 


0.1600 


0.09000000000000 


0.05500000000000 


0.03750000000000 


NC 3 (x s ) 


0.3000 


0.1600 


0.10200000000000 


0.07511111111111 


0.06111111111111 


Go (x s ) 


0.3000 


0.0900 


0.02700000000000 


0.00810000000000 


0.00243000000000 


Gi(x s ) 


0.3000 


0.1600 


0.10200000000000 


0.06725714285714 


0.04459836734694 


Wi(x s ) 


0.3000 


0.1600 


0.09000000000000 


0.05080000000000 


0.02868000000000 



Notice that this error estimate in the Lebesgue case is pessimistic because 

we have L f(x)dx — ^(f) = — for suitable £. We will obtain the 

2880 

optimal constant in the error estimate in the next remark when considering 
the fact that this rule coincides with an Hermite rule. 

Notice also that this quantity is always null only in the cases a = 0, 1,4 that 
gives exactly the information that the Newton- Cotes three points-formulae 
is an exact formulae only in the dirac cases. □ 

Remark 3.3 (Error Estimate 2: GL2). The procedure used above for the 
estimate of the error is not useful in the case of GL2 because we know, by 
construction, that = K~ , but it is not easy to compute this value, due 
to the fact that it is given by an integral on a non dyadic subinterval. We 
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Table 5: Calculated moments, case a = 0.45 





s= 1 


o 

Z 




4 


r 

5 


Jo •£ "^0.45 


4^00 

U.4:UUU 




90^0^71 49&^71 


1 fi00949&^71 49Q 


1 ^0071 4fi^1 3**fU 

\J. ±0\J\J t 14U0100U1 


NC (x s ) 


0.5000 


0.2500 


0.12500000000000 


0.06250000000000 


0.03125000000000 


NCi(x s ) 


0.4500 


0.4500 


0.45000000000000 


0.45000000000000 


0.45000000000000 


NC 2 (x s ) 


0.4500 


0.2850 


0.20250000000000 


0.16125000000000 


0.14062500000000 


NC 3 (x s ) 


0.4500 


0.2850 


0.20603571428571 


0.16373809523810 


0.13898809523810 


(Go(x s ) 


0.4500 


0.2025 


0.09112500000000 


0.04100625000000 


0.01845281250000 


&i(x s ) 


0.4500 


0.2850 


0.20603571428571 


0.15456581632653 


0.11703565233236 


Wi(x s ) 


0.4500 


0.2850 


0.20250000000000 


0.14805000000000 


0.10894500000000 



can obtain an error estimate, although, noticing that the integration rule 
corresponds to integrate the Hermite polynomial with the central node of 
multiplicity 2. Indeed, see f2l\j %8.4, the interpolating polynomial in this 
case can be written as: 

rr (\ ^f(v\\ i \ . ( x ~ Co)(^ ~ Ci)(^ ~ C2) ^ 

U f,H{x) =2^f{Q\ q {x) + — — — / (Cl) 

(Ci — C0KC1 — C2J 

where the functions X q (x) are the usual Lagrange polynomials. Notice that 
the function that multiplies the /'(Ci) is the nodal polynomial, and this im- 
plies that the integral of this term gives no contribute. By the other side, we 
can apply the usual interpolation error estimate for the Hermite interpola- 
tion, that in our case gives: 

finite \ 

fix) - n />H (x) = (x - Co)0 - (i) 2 (x - C2) for some & € [0, 1] , 

see \2tf %8.4- This leads, for the error in the GL2 case, to: 
-1 ri fimte \ 

f(x)dn a -GU(f)= / ^P^-l)(x-Cl) 2 ^ a . 

10 Jo 4! 

From this, explicit calculation gives: 

2a(17a 3 - 34a 2 + 9a + 8) /""(fl 

24— f° r some Z e [0, 1] , 

where we have applied proposition ^. 1\ to the function f""(£ x ) as ^ n (|6fa| . □ 
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3.2 Composite rules 



The usual way in which the quadrature rules are used in composite manner 
is to introduce a partition of the initial interval and to consider on each 
subinterval the integral to be approximated with a proper quadrature rule. 
In this paragraph we will see how to do this in the framework of integration 
with respect to binomial measures. 

We introduce a partition of the initial interval in N subintervals Ji = 
such that = xq < x\ < • • • < xjy = 1. We know how to rescale 
integrals on dyadic intervals -by means of lemma 12.51 and for this reason we 
consider the next definition. 

Definition 3.4 (Dyadic-regular Partitions). We will say that the partition 
{Ji}i=o...N~i is dyadic-regular if\/i3j* and k* s.t. Ji = Xj„ . 
We will say, in particular, that it is dyadic-proper if N = 2 k and Ji = 
X? Vi = 0,...,2 fc -1. 

On each subinterval we consider the function to be approximated with a 
p-polynomial. As seen in the previous chapter, we can use different choices 
of interpolating polynomials, leading to different quadrature rules; we will 
consider that in each subinterval we apply the same interpolatory rule. The 
corresponding quadrature formula on the subinterval will be called local 
quadrature rule. In our notation ]L (/) will indicate that we are applying 
the local quadrature rule I p on N subintervals to the function / £ Lz a : 

N-l V 
i=0 q=0 

When we will write Ip (/) we will consider the dyadic-proper case. 
In order to write local quadrature rules we modify nodes and weights seen 
in tables [T][2j In particular it is easily seen, by lemma 12.51 that if the the 
p + 1 nodes Q , q = 1 . . . p + 1 on Jj = Xj» are taken simply by rescaling 
them from [0, 1] in the interval Ji (C = (j* + Cg)/(2 fc )), the corresponding 
weights are /3* = j3 q ^i a {Ji). When we apply the formula only on one dyadic 
interval we will use the notation Ji). 



Remark 3.5 (Extrapolation). We have already noticed that both the NC2 
and GL2 reduce to the 3 points Newton- Cotes formula in the Lebesgue case. 
In this case (a = 1/2) we have that the 5 points formula can be written as 
(see m$): 

NQ/ = 16(NC 2 {f, Xl ) + NC 2 U,X\))- NC 2 (f, X° ) E 

15 1 

16(GL 2 (/, Xl) + GL 2 (/, X, 1 )) - GL 2 (/, gg) E 
— ^ J 
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For general 1/2 the rule NCf uses 5 points and is of degree of exactness 
2, while the rule GL^ uses 6 points and is of degree of exactness 4, thus none 
of the two are of interpolation type. On the 6 points of the rule GLf' we can 
construct the rule of interpolation type that is of order 5. This formula also 
coincides with the Newton-Cotes formula NC4 in the Lebesgue case. We will 
denote this formula with M4. Nodes and weights are: 



c 

A) 
Pi 
Pi 

Ps 
Pa 
Ps 



0. 



3a + 2 



1 1 

2'2 + 



3a + 2 , 3a + 2 



1, 



14 '2' 2 14 
688a 5 - 24257a 4 + 59238a 3 



7 

32825a 2 



19584a + 16740 



1395(a + 3)(3a + 2) 2 
224a(1667a 4 - 6012a 3 + 4855a 2 + 1442a - 1952) 

1395(4 - a) (3a - 5) (9a 2 + 12a + 4) 
32a(43a 3 - 86a 2 - 669a + 712) 



1395(3a 
224a (1667a 4 - 



-2)(5- 
2323a 3 



3a) 



2523a 2 + 3395a - 216) 



1395(3a - 5) (9a 3 + 18a 2 
a(688a 4 + 20817a 3 - 30910a 2 



37a - 30) 
6227a - 1108) 



1395(a-4)(3a-5) 2 
98a(2311a 3 - 4622a 2 + 1137a + 1174) 

162a 3 



1395(81a 4 - 162a 3 - 99a 2 + 180a + 100) 

Note that the first five nodes are in increasing order, while the last £5 co- 
incides with the midpoint in the Lebesgue case and is between d and £2 if 
a < 1/2 and between £2 and £3 in the other case. □ 

Note that the estimate (16flp can be written on a dyadic subinterval by 
the following: 



f f(x)dfi a -I p (f,X*) 
Jx* 



M4) 



2*(p+i)(p+ 1)! 



vfc yk vk 



co p (x) djj, a 



(9) 



By ([9]), taking into account that an analogous equation is true for K a , we 
obtain, for suitable 771,772 £ Xf: 



[ f(x)dfi a -I p (f,X 
Jx* 



< 



< 



2*(P+l)(p+ 1) 



min{if+ ,K a } 

+ f {p+1) (m 



f ip+2) (vi) 
1 



+ 



u) p {x) dfj, c 



(Hb) 
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Moreover, we can write in the single dyadic interval also the error estimate 



f(x)dfx a -I(f(x),Xf] 



(r + 1)! 

(r + 1)! 



r+1 



iV 



^) dfi a -Y J ^(x^)p q 

q=0 



Cq-j 
2 k 



J_ 

2 k 



r+1 



A" 



g=0 



r+1 
9 



r+1 



(10) 



The most simple way to consider the composite integration rules is to 
consider dyadic proper partitions of increasing order. We report a simple 
algorithm that allows us to introduce the definition of order of convergence. 
We consider the following iterative procedure: 

Composite Algorithm 



Initialization 
Cycle Control 
Local Quadrature Application 
Bisection 



Output: 



put k = 

while STOP CRITERION 

Compute l(f, X k ) = lf ] Vj = 1 . . . 2 k 

k = k + l 
end while 

_ ^2 fc T (fc) 



Routine utilized: Local quadrature rule I 



We have written the scheme thinking at a procedure that calculates an error 
estimate and iterates until the result is considered satisfactory: the stopping 
criterion has to be chosen to complete the scheme. Usually it is considered 
in order to satisfy an error requirement: 



fdfi c 



< tol 



It can also be chosen to avoid too many numerical computations. 

For this algorithm we are interested in convergence properties. We will 
say that the composite rule converges of order 7 if 



f(x)d» a -I 2 p (/) 



for / sufficiently regular. Applying the Taylor expansion of the function with 
Peano's remainder ()10p we obtain that the formula has order of convergence 
at least equal to the degree of exactness. 
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Proposition 3.6. Let I be of degree of exactness r and consider f to be 

r + 1 times derivable. 

Then I has order of convergence at least r. 

Proof: From the error estimate (I10D we have: 



/ f(x)d f i a -I(f,Xf) 
Jx* 



(r + 1)! 



r+1 



„r+l 



dfl c 



N 
q=0 



r+1 

q 



Now, summing all the subintervals: 



f(x)d^ a -lfU) 



7 (r+1) (0 

(r + 1)! 



k ar t , for some £ £ [0, 1] (11) 



where we have called k a)r j = Jq x r+1 dfi a — Y1^=q PqCq +1 ■ This estimate 
gives the requested property. □ 
In figure © we plot -for the same two test function seen in the previous 
section and in the case a = 0.3- the convergence of the composite algorithm. 
From this plot we can see that the order prescribed by the previous theorem 
is confirmed. 

Remark 3.7 (Stopping criterion for the non adaptive algorithm). This con- 
vergence property can be used also to write a stopping criterion based on an 
error estimate, when it is available an estimate of the p-th derivative of 
the function. The most common strategy, although not very reliable (see 
13 HE/), is to estimate the maximum of the derivative considering the in- 
formation known from equation (|11|) at level k and k + 1. Starting from 
relation (|10p on a single interval Xj we can write that: 



L 



fdn a -I(J,X} 



r+1 



M*ft/ (r+1) (e*' fc ) 

(r + 1)! 



and 



f fdna-W'X&^ + Hf,*: 



k+l 
2j+l 



)] 



+ 1 MQ (X 2 Y 1 )/^ +1 )(^V^(^+\)/ (r+1 H^ 1 ) 



fc+i 



2 k+l 



(r + 1)! 



ka,r,] 



If we consider /( r+1 ) almost constant in X*- to the value Kj'l 1 
sider to approximate this value with: 

2 (fc+l)(r+l)( r + 1 

(2-+ 1 - l)fi a (Xf)k 



we can con- 



K r+1 



9(i+l)(r+l)/.ii\| ( 
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Figure 3: Plot of the convergence history for some of the integration rules. 
Test functions are in equation (JS|) . 



Define, now, K'L +1 



L fc - max j=0,...,2 K -l 

consider as stopping criterion in the non- adaptive composite algorithm: 



K r+l 



With these positions we can 



kmin < k < min 



1 kgjKj^l 
r ° g2 (r + l)\tol 



+ 1, fc r , 



(12) 



where [•]+ denotes the integer part. 

Notice that k m i n is considered to force the algorithm to do the first com- 
putations (ad example in the case of peaked functions) and k max to avoid 
too many computations. 



□ 



4 Final Remarks 

We conclude with some remarks and conclusions. 

We develop quadrature formulae for the family of binomial measures. More- 
over error estimates, based both on interpolation errors and on Taylor ex- 
pansion with Peano's remainder, for such formulae have been established. 
We eventually list some possible future developments. 

1. Develop an adaptive algorithm. 
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2. Consider how to calculate numerically rules of higher order. 

3. To extend the analysis to more general classes of measures, with special 
attention to measures given from experimental data. 
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